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A theory for (l+3)-dimensional relativistic Brownian motion under the influence of external force 
fields is put forward. Starting out from a set of relativistically covariant, but multiplicative Langevin 
equations we describe the relativistic stochastic dynamics of a forced Brownian particle. The cor- 
responding Fokker-Planck equations are studied in the laboratory frame coordinates. In particular, 
the stochastic integration prescription, i.e. the discretization rule dilemma, is elucidated (pre-point 
C**) ' discretization rule vs. mid-point discretization rule vs. post-point discretization rule). Remarkably, 

within our relativistic scheme we find that the post-point rule (or the transport form) yields the only 
Fokker-Planck dynamics from which the relativistic Maxwell-Boltzmann statistics is recovered as 
the stationary solution. The relativistic velocity effects become distinctly more pronounced by going 
from one to three spatial dimensions. Moreover, we present numerical results for the asymptotic 
. mean square displacement of a free relativistic Brownian particle moving in (1+3) dimensions. 

PACS numbers: 02.50.Ey, 05.40.-a, 05.40.Jc, 47.75. +f 
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The problem of formulating a consistent theory of Brownian motions 0, 0, ES 0, IE in the framework of special 
relativity 0, 0] represents a longstanding issue in mathematical and statistical physics (classical references a re Iq |. fiol 
lllj : more recent contributions include (30,0,0,0,0,0,0,113; for a kinetic theory approach, see |2ltl22j1. 
In a preceding paper |23j | - referred to as paper I hereafter - we have discussed in detail how one can construct 
Langevin equations for (l+l)-dimensional relativistic Brownian motions. In particular, it was demonstrated that 
for the relativistic Langevin equation per se cannot uniquely determine the corresponding Fokker-Planck equation 
(FPE). This ambiguity arises due to the fact that the relativistic Langevin equations, when e.g. written in laboratory 
£SJ ' coordinates, exhibit a multiplicative coupling between a function of the momentum coordinate and a Gaussian white 
C*") \ noise process (laboratory frame = rest frame of the heat bath). Thus, depending on the choice of the discretization 
rule, different forms of relativistic FPE are obtained [3, Ell- 
in paper I, we have analyzed the three most popular discretization rules for Langevin equations with multiplicative 
noise, which can be traced back to the proposals pioneered by Ito (pre-point discretization rule prl |2?| ) , by Fisk and 
Stratonovich (mid-point rule [25|, |28|, |2!J H3, H3 ) , and by Hanggi and Klimontovich (post-point rule |32l l33t l34i l35| ) . 
In this context it could be shown that only the Hanggi-Klimontovich (HK) interpretation of the Langevin equation 
yields the transport form of the Fokker-Planck equation with state-dependent diffusion, whose stationary solution 
coincides with the one-dimensional relativistic Maxwell distribution. The latter is known from Jiittner's early work 
I ■ on the relativistic gas [3(| |37j and also from the relativistic kinetic theory • 

In paper I, we have focussed exclusively on the simplest situation, corresponding to free (l+l)-dimensional rel- 
ativistic Brownian motions. Therefore, the present paper aims to extend the analysis to the physically relevant 
(l+3)-dimensional case. In particular, we wish to include as well the effects of additional, external force fields. To 
this end the paper is structured as follows: In Sec. [H] the relativistic Langevin equations are given in covariant 
4-vector notation and also in laboratory frame coordinates. The corresponding Fokker-Planck equations and their 
stationary solutions are considered in Sec. II I II Section llVI contains a discussion of numerical results for the mean 
square displacement of free Brownian particles. The paper concludes with a resume of the main results in Sec. 
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II. RELATIVISTIC LANGEVIN DYNAMICS 



Let us first discuss the manifestly Lorentz-covariant 4- vector form of the relativistic Langevin equations (Sec. Ill Afl . 
For that purpose, we shall use of the results derived in Sec. II of paper I, which can readily be generalized to (1+3)- 
dimensions. Subsequently, the relativistic Langevin equations will be written in laboratory coordinates (Sec. IIIBj l. 
The latter form provides the basis for the numerical results of Sec. IIVI 

With regard to notation, the following conventions will be used throughout the paper: Upper and lower Greek 
indices a,/?,... can take values 0,1,2,3, where '0' refers to the time component. Upper and lower Latin indices 
. . . take values 1,2,3 and are used to label the components of spatial 3-vectors, denoted by bold symbols. For 
example, we write (x a ) — (x°,x) = (ct,x l ) and (p a ) — (p°,p) = (E/c,p l ) with t denoting the coordinate time, E 
the energy, c the vacuum speed of light, x l and p l the spatial coordinates and relativistic momenta, respectively. 
Moreover, Einstein's summation convention is applied throughout. The (l+3)-dimcnsional Minkowski metric tensor 
with respect to Cartesian coordinates is taken as 

iVap) = (V af3 ) = diag(-l, 1, 1, 1) , ( V a p ) = ( V J) = diag(l, 1, 1, 1). 

As commonly known, covariant vector components x a can be calculated from the contravariant components x a by 
virtue of x a — rj a px^ , which, in particular, means that for Cartesian coordinates xo — —x° and Xi = x l holds. 
Further, if in a certain inertial coordinate system S a Brownian particle has the 3- velocity v(t) = dx(t)/dt, then the 
differential dr of its proper-time is defined by 

dr = d^l-^. (1) 



A. Langevin equations in 4-vector notation 



Consider a Brownian particle with rest mass m, proper time r and 4- velocity u^(t), i.e., the 4- momentum of 



the particle is given by p a = mu a , where u a u a 



Assume that the particle is surrounded by an isotropic, 



homogeneous heat bath with constant 4- velocity U 13 and, additionally, subject to an external 4-force K a [x v , p^) such 
as, e.g., the Lorentz force. Then, according to the results in Sec. II of paper I, the relativistic Langevin equations of 
motion read 



dx a [r) = — dr 
m 

dp a ( T ) = {K a -v a [p -mU p ]} dr + w a ( T ). 
For an isotropic homogeneous heat bath, the friction tensor v a p in i|2bj) is given by 

u a up " 



v"0 = v\ V"p + 



(2a) 
(2b) 

(2c) 



with v denoting the scalar viscous friction coefficient measured in the rest frame of the particle. Furthermore, the 
relativistic Wiener increments w a (t) = dW a (r) are distributed according to the probability density 



V i+d [w a (T)} = 



{AttD dr) 3 / 2 



exp 



w a (r) w a (r) 
4Ddr 



X S\ 



(2d) 



where D is the scalar noise amplitude parameter measured in the rest frame of the particle. Some useful comments 
concerning Eq. (J2J are appropriate: 

(i) The covariant friction tensor in Eq. Ij2c|l carries the same structure as the covariant pressure tensor for an 
ideal fluid. In particular, this means that in each instantaneous rest frame £* of the particle, where temporarily 
(m™) = (c, 0) holds, the tensor form Eq. lj2"c)) reduces to the diagonal form 



/0 ON 
Of 00 
v 

\0 v) 



(3) 



This form of the friction tensor reflects the simplifying assumptions that the heat bath can, in good approximation, 
be considered as isotropic and homogeneous. 
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(ii) Note that the probability density of the increments in Eq. I|2dj) can equivalently be written as 



V 1+3 \w a ] = 



(4ttL> dr) 3 / 2 



exp 



1 

2d7 



x 5 [u a w a ] , 



where the tensor 



D f>- 2D 



(4a) 



(4b) 



carries the same isotropic structure as the friction tensor from Eq. lj2*c| . The 8- function in Eqs. (|4a|l accounts for the 
fact that the Minkowski scalar-product of the 4-force and the 4-velocity must identically vanish, 
(iii) The increment density in Eq. (|2d|l is normalized so that 



(5a) 



{ 3 p oo | 
J] / dK«] \V 1+3 [w a (r) 



holds. For the first two moments one finds 



K(r)) = 0, 



0, r ^ r'; 

D^dr, t = t\ 



where 



(5b) 
(5c) 

(5d) 



The easiest way to validate this is to perform the calculations leading to Eqs. JSJ in a co-moving Lorentz frame E», 
where, at a given instant of time i»(r), the particle is at rest. In such a co-moving frame £* the marginal distribution 
of the spatial momentum increments, defined by 



-oo 



(6) 



reduces to a Gaussian. One thus recovers, as it should hold true, from Eqs. (J2J the nonrelativistic Brownian motion 
in the Newtonian limit case v 2 <§; c 2 . 



B. Langevin dynamics in the laboratory frame 

In this section the covariant Langevin equations J5J will be rewritten in laboratory coordinates. A laboratory frame 
£o is, by definition, an inertial system, in which the heat bath is at rest. That is, in So the 4-velocity of the heat 
bath is given by (U 13 ) — (c, 0) for all times t, where t is the coordinate time of So- 

From Eq. (|2a|> . we obtain three differential equations for the position coordinates 

dx i (t) = v i dt, (7a) 

where 



(7b) 



^Jm 2 c 2 + Pip 1 

Furthermore, since we have (U 13 ) = (c, 0) in So, the four stochastic differential Eqs. i|2bjl can be rewritten as [23l |3S| 

dp 1 = {~t~ x K l -up l )dt + w\ (7c) 
dE = (j^K 1 -vp l )v t dt + cw°, (7d) 

where the relativistic (kinetic) energy is here defined by E = cp°, and 

/ v-v { \ ~ 1/2 ( p ■ i \ 1/2 
^ V c 2 J \ m 2 c 2 J ^ ^ 
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Before discussing the stochastic increments w a , let us briefly consider the deterministic force components K l . If 
F = (F l ) is the nonrelativistic (Newtonian) force, then the corresponding relativistic force 3- vector K — (K l ) is 
given by (see Chap. 2 of |38j ) 



(v ■ F) 

K = F + { 1 ~l)^ r ±v. 



(8) 



If the dynamics is confined to one spatial dimension only, then Eq. © simplifies to K = jF. 

Next, let us turn to the stochastic force components, appearing on the rhs. of Eqs. (JTcJ and (|7d|l . According to 
Eq. il2d(l . the distribution of the stochastic momentum increments w a , also depends on the particle's velocity v. By 
virtue of the relation 



(u a ) = (-7c, jVi) , (w a ) = (w°,w l ) 

we can rewrite the increment density (I2d|l in laboratory coordinates as follows 



P 1+3 \w a ] 



I 7 \ 3 /2 

(^Ddt) 6XP 



AD dt/j 



x S (c7W° — jViW 1 ) 



(9) 



(10) 



As we already pointed out above, the (5-function in Eq. (|10|l reflects the fact that the energy increment w° is coupled 
to the spatial momentum increments w l via 



= u a w a = —cjw + jViW 1 => w = 

This is just the stochastic analogue of the well-known deterministic identity 

u a K a = 0. 

Hence, similar to K° , also w° can be eliminated from the Langevin equation |J2 

dE = (j^K 1 ~vp l )v t dt + v lW l {t) = v % dp 
Using the identity (|7b|l . we can further rewrite Eq. I|12|l as 

CPi 



(11) 



yielding 



dE 



■dp 1 , 



^/m 2 c 2 + p 2 c 2 

where p 2 = Pip 1 ■ From the preceding equation we regain the well-known energy momentum law 



E(t) = v/m 2 c 4 + p(t) 2 c 2 



(12) 



(13) 



(14) 



It might be worthwhile to remark that in the presence of force fields the relativistic energy E = cp° is generally 
different from the relativistic Hamiltonian [3|j . Note that Eq. i|14|) remains valid also for our stochastic model. 

The relativistic Brownian motion is therefore completely described by the three Langevin equations (|7c|l . The 
statistics of the increments w l in l|7cf) is determined by the marginal distribution V [w], defined in Eq. (jSJ. Performing 
the integration over the 5-function in Eq. (|10fl , we find 



V 3 [w] 



dw° P 1+3 \w a ] 



— OO 

oo 



dw° 
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\4TrDdtJ 



AttD dt 

3/2 



3/2 



exp 



exp 



7 



ADdt 



4D dt/j 

ViVi 

7 2 



8 (c7w° 



(15) 



where dij denotes the Kronecker-symbol (defined by Sij = 1 if i = j, and 5y = otherwise). 

In principle, it is straightforward to perform computer simulations on the basis of Eqs. lf7cfl and (|15|) . In Sec. IIVI 
we will discuss several numerical findings. Before doing so, however, it is worthwhile to consider in greater detail the 
Fokker-Planck equations of the relativistic Brownian motion in the laboratory frame E n. B y d oing so it will become 
clear that the choice v l = v l {t) in Eq. i|15|) 1S consistent with an Ito-interpretation |24l 12a. |27| of the stochastic 
differential equations l|7c|) . However, we shall also see that alternative interpretations - such as e.g. the post-point 
discretization rule, v % = v l (t + dt) - lead to physically reliable results as well. 
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III. RELATIVISTIC FOKKER-PLANCK EQUATIONS 



The objective here is to discuss relativistic Fokker-Planck equations (FPE) for the one-particle momentum density 
f{t, p) and, as well, for the phase space density f(t,p, x). In the remainder, we will exclusively refer to the coordinates 
of the laboratory frame So- Before turning to the relativistic FPE in Sects. IIII Bl and lHI CI it is useful to briefly recall 
the nonrelativistic case. 



A. Nonrelativistic case 



Consider the nonrelativistic Langevin equations 



dx l = v l dt 

dp 1 = [K r -vp l )dt + w l 



(16a) 
(16b) 



where p l (t) — mv l {t) denotes the nonrelativistic momentum components, K % = —d l U represents the vector com- 
ponents of a conservative force with time-independent potential U(x), and where the increments w 1 = dW 1 are 
distributed according to 



V[w] = 



1 



iirDdt 



3/2 



exp 



Wi w 



AD dt 



(16c) 



Equations l|16|) govern the nonrelativistic motions of a Brownian particle in the rest frame of the heat bath. As one 
readily observes, in the case of conservative force fields, Eqs. (|16|) can be obtained from the relativistic equations 
(|7a|l , l(7c|) and (|15|) by formally taking the limit c — > oo (Newtonian limit case) . It is well known that the phase space 
density f(t,p,x), associated with the stochastic process (|16|) . is governed by the FPE p3. 125L liol] 



Ot 771 ox op 1 op 1 \ op 1 

The stationary solution of Eq. I|17|l is the (nonrelativistic) Maxwell-Boltzmann distribution, i.e., 

p 2 + 2m U(x) 



f{p,x) = Cexp 



2m k^T 



(17) 



(18a) 



where C is a normalization constant. The temperature T of the bath is defined by the Einstein relation (fcg denotes 
the Boltzmann constant) 



k*T = 



D 



The related marginal momentum distribution is the usual Maxwellian probability density 



2ttD 



3/2 



exp 



up- 
2lJ 



(18b) 



(19) 



B. Relativistic FPE for free Brownian particles 



We next reason three different types of relativistic Fokker-Planck equations for the momentum density f(t,p) of 
a free Brownian particle, whose dynamics is governed by the stochastic process given by Eqs. (|7cf) and i|15|) . The 
corresponding equations for the phase-space density f(t,p,x) will be considered separately in Sec. IIII CI where we 
will also include the influence of an external force field. 

As our starting point serves the relativistic Langevin equations JJcJ), which hold in the laboratory frame So (i.e., 
in the rest frame of the heat bath) . In absence of external force fields the stochastic differential equations (I7c)l reduce 
to 



dp 1 = -h>p l dt + w\ 



(20a) 
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where p 1 = ^mv 1 is now the relativistic momentum and 7 has been defined in Eq. l(7e|l . According to Eq. I|15|) - and 
in distinct contrast to Eq. (|16cl) - the distribution of the relativistic increments reads 



7 KAnDdt, 
with matrix elements given by 



v3/2 




J exp 


ADdt 



(20b) 



Following the reasoning of paper I, the next aim is to rewrite the Langevin equations (1201) in such a form, that the 
resulting equations exhibit multiplicative Gaussian white noise, governed by a velocity-independent normal distribu- 
tion of the form (|16cl) . In order to achieve this objective we first note that the matrix A(p) = (A l j) is symmetric. Its 
eigenvalues and determinant are given by 

spec(A) = {7, 7, 7' 1 }, det(A) = 7. (21) 

Thus, the matrix A is positive definite for velocities v 2 < c 2 , and the elements of the inverse matrix A^ 1 read 



'Sh | v 3 Vk \ _ ( ci , P 3 Pk \ 1 
7 2 c 2 J 1 

Furthermore, there exists a unique Cholesky-decomposition 



k \ 7 2 C ? ) \ m 2 c 2 J 7 




,L 3 i L 3 2 L 3 3/ 

where the matrix L{p) is nonsingular with elements given by 



<L\ L 2 , L^i 

A = L 1 L = I L 2 i L 2 2 II L 2 2 L 3 2 I , (23) 



(24) 





= VAh, 


L\ 


= A\/L\, 


L 2 2 


= VAh- 


L\ 


= A\/L\, 


L\ 


= {A\-L 


L\ 


= VAh- 



3 T 2 



The inverse matrix L(p) 1 reads 



L 2 2L 3 3 —L 2 iL 3 3 L 2 \L 3 2 — L^iL 2 ^ 
L 1 = I L3 3 £l i ~L\L\ ) , (25a) 



det ( L ) I " IT 1 L\LK 



where 

det(i) = L\L 2 2 L 3 3 . (25b) 
Let us next introduce a stochastic vector variable y(t) = (y l (t)) by 

//' /.',«•'. (26) 
Then, by taking into account that for u; T = (wi) and y T = (?/i) the relation 

w T Aw = w T L T Lw = (Lw) T Lw = y T y (27) 
holds, we can rewrite the Langevin equation l|2U|) in the form 

dp = -vpdt + L(p)~ 1 y, (28a) 
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where y(t) is distributed according to the momentum- independent Gaussian probability density 

Vi 



4L> dt 



(28b) 



Put differently, because the inverse matrix L(p)^ 1 depends on the momentum coordinate p, the random vector y(t) 
enters in relativistic Langevin equation (|28a[l as an ordinary 'multiplicative' Gaussian white noise process with noise 
strength D. 

As is well known, for multiplicative stochastic processes of the type l|28fl the Langevin equation per se does not 
uniquely determine a corresponding Fokker-Planck equation 0>E3- In the following subsections, we shall discuss the 
three most popular choices of resulting Fokker-Planck equations for a Langevin equation of the form 128|). These choices 
are rooted in the different proposals put forward by Ito [3 HE HE H3 > by Stratonovich and Fisk |25l EsL I2H EfflL l3ll| . 
and by Hanggi [5E HE 13 an d Klimontovich j^E , respectively. All three approaches have in common that the related 
Fokker-Planck equation can be written as a continuity equation (conservation of probability) of the form [34[ 

|/(* ) p) + | I i i (i,p)=0, (29) 

but with distinct different expressions for the probability current j{t,p). It is worthwhile to anticipate here that only 
the Hanggi-Klimontovich approach (see Sec. IIII B 3|) yields a stationary distribution, which can be identified with 
Jiittner's relativistic Maxwell distribution |36j . 

1. Ito approach 

According to Ito's interpretation of the Langevin equation (|28a|l . the coefficient matrix before y(t) is to be evaluated 
at the lower boundary of the interval [t, t + dt], i.e., 

L(p)- l = L{p{t))-\ (30) 

Ito's choice is also known as the pre-point discretization rule 0> HE HE H3 an d leads to the following explicit 
expression for the current 



(31) 



In view of the identity 



L- 1 {L- 1 ) T = L _1 (L T ) _1 = [L^L)- 1 = A~\ (32) 
Eq. (|31|l can be rewritten more conveniently as 

d 



yp'f + D^-iA^fJ 



(33) 



where the matrix A(p)~ 1 is given in Eq. The related relativistic Fokker-Planck equation is obtained by inserting 

this current into the conservation law (|29|l . As elucidated in the Appendix lAl the stationary solution of the resulting 
Fokker-Planck equation reads 



p 2 v 3/2 ( I p 2 



/i(p) = Ci 1 + ^0 cx p + , (34) 

V m z c J V V m z c I 

which contains in the prefactor an explicit (non-Maxwell like) dependence on velocity and where C\ is a normalization 
constant. In the solution (|34|l . the dimcnsionless parameter 

X = D (35) 
is related to the scalar temperature T of the heat bath via the Einstein relation 

, „ mc 2 D 

k B T= = . 36 

X mv 

Thus, the quantity \ = mc 2 / (k^T) gives the ratio between rest energy and thermal energy of the Brownian particle. 
It should be mentioned that we used in paper I the notation (3 instead of x- However, in order to avoid a possible 
confusion with the commonly used abbreviation /3 — {k^T)~ l we opted here for this slight change of notation. 
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2. Stratonovich-Fisk approach 

According to the Stratonovich-Fisk prescription, the coefficient matrix before y(t) in (|28a|l is to be evaluated with 
the mid-point discretization rule, i.e., 

L(p)-i=L(M + « t + dt A~\ (37) 



V 2 

This choice [HHHHH^] leadS t0 & different ex P ression for thc current pi El l30| : 



The explicit stationary solution of Stratonovich's Fokker-Planck equation reads (see Appendix IA 21 



P 2 \ '' ' ( P 2 



f B *<P) = C BE \l + ^) exp^l + ^J, (39) 
and thus differs from Eq. (|34|l in the power of the velocity-dependent prefactor. 

3. Hanggi-Klimontovich approach 

Ultimately, let us next consider the Hanggi-Klimontovich (HK) stochastic integral interpretation, sometimes referred 
to as the transport form [32l l33l l34 | or also as the kinetic form (35(. According to this interpretation, the coefficient 
matrix in front of y(t) in (|28a|l is to be evaluated at the upper boundary of the interval [t,t + dt]; i.e., within the 
post-point discretization we set 

Lip)- 1 =L(p(t + dt))-\ (40) 
This choice leads to the following expression for the current (3^, 03 HU 



Jhk(*>P) 



3 dp j . 



(41) 



and the stationary solution of the related FPE reads (see Appendix Eq. (|A 3|) ~) 



/hk(p) = C* HK cxp \ - X \j 1 + ■ ( 42a ) 

Note that this solution contains no velocity dependence in the prefactor. Using the temperature definition in 136|) 
and the relativistic kinetic energy formula E = m 2 c A + p 2 c 2 , we can recast (|42al) as 

/hk(p) - C HK exp (--^=) . (42b) 



. k B T 

The normalization constant is given by [TlT | 



C H k= / d 3 P cxp( - X \}l + J!L ) =4^ / dpp 2 exp(- X Jl + -^) = 47r(mc) 3 ^^, (42c) 



where K 2 (x) denotes the modified Bessel function of the second kind. The distribution function (|42|l is known as 
the relativistic Maxwell distribution. It was first derived by F. Juttner [3?3 | in 1911, when he investigated the velocity 
distribution of non- interacting relativistic gas particles (see also [37J ). By comparing (I34|l . (|39f) and l|42a|) one readily 
observes that the stationary solutions /i/sf differ from the Juttner function /hk through additional p-dependent 
prefactors. The quantitative difference between these three stationary solutions becomes significant in the relativistic 
limit, corresponding to a low rest energy-to-temperature ratio \- 
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As already mentioned in paper I, for the one-dimensional case the relativistic Maxwell distribution has also been 
obtained by Schay, see Eq. (3.63) and (3.64) in Ref. ||, who studied relativistic diffusions employing a transfer 
probability method. Moreover, the distribution l|42|l can also be derived in the framework of the relativistic kinetic 
theory [22]. This suggests that the HK-discretization rule is physically preferable, if one wishes to use the above 
Langevin equations for numerical simulations of relativistic kinetic processes. In general, however, additional infor- 
mation about the microscopic structure of the heat bath is required, in order to decide which discretization rule is 
physically reasonable (see, e.g., the discussion of Ito-Stratonovich dilemma in the context of 'internal/external' noise 
as given in Chap. IX. 5 of van Kampen's textbook |24|). 

Finally, we still note that the related velocity probability density functions 0i/sf/hk('w) are obtained by applying 
the general transformation law 



(«) = f(p(v)) 



dp 



dv 



(43) 



where, as usual, 



p(v) = 



The determinant factor 



dp 



dv 



TYl 1 



2\ -5/2 



(44) 



appearing on the rhs. of Eq. I|43|) . ensures that the velocity density functions ^i/sf/hk(' 1 ') drop to zero for v* 
For example, in the case of the Juttner function (|42|l one explicitly obtains 



^hk(«) 



2\ -5/2 



4irc 3 K 2 (x) 



exp 



y/1 - V 2 /C 2 



(45) 



C. The inclusion of external force fields 



The preceding section concentrated on Fokker-Planck equations for the momentum density f(t,p). In this part we 
shall discuss the corresponding resulting equations for the one-particle phase space density f(t,p,x). As before, we 
refer to the coordinates of the laboratory frame So, in which the heat bath is at rest. 

If an external force field is present, then Eq. (|28a|l generalizes to 



dp = {^K - vp)dt + L(pY 



(46) 



where y is distributed according to the momentum-independent Gaussian density from Eq. 
relativistic Fokker-Planck equation for the full phase space density f(t,p, x) thus reads 



The related 



l f+ pL± f + ±(^ f 

dt mj dx l dp 1 \ 7 



d_ 

dp 



T^I/SF/HK 



(47) 



where the current densities Ji/sf/hk(*'P' x ) are obtained by replacing f(t,p) with f(t,p,x) in the above expressions 
for Ji/sf/hk(*iP)j respectively. In the limiting case that v — > 0, D — > 0, the rhs. of Eq. (|47|l becomes equal to 
zero and one regains the relativistic Liousville equation or, equivalently, the collision-less Boltzmann-Vlasov equation 

E1E2. 

As a particular example, let us consider a relativistic Brownian particle with rest charge q, being subject to a static 
electromagnetic field (E,B), measured with respect to So- Then, in addition to the stochastic interaction with the 
heat bath, the deterministic Lorentz-force Hal 



K(p,x) 



E(x) 



P 



x B(x) 



(48) 
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is acting on the particle, where 'x' denotes the exterior vector product. For simplicity, let us confine ourselves to the 
HK-form of the FPE and let E(x) = V$(a;) and B(x) = in the laboratory system. In this case, the stationary 
solution of Eq. lj4*?|) emerges as 



fuK(p,x) = C H Kexp 



= C H Kexp 



lip) ~ 



E(p) - q$(x) 



(49a) 
(49b) 



where Chk is a normalization constant, and E(p) = (m 2 c 4 + p 2 c 2 ) 1 ^ 2 denotes the relativistic (kinetic) energy. It is 
reassuring to see that the solution (|49|l just represents the relativistic generalization of the nonrelativistic Maxwell- 
Boltzmann distribution from Eq. 



IV. NUMERICAL INVESTIGATIONS 



The numerical results presented in this section were obtained on the basis of the relativistic Langevin equations J2SJ|, 
which hold in the laboratory frame So. For simplicity, we confined ourselves to considering free Brownian particles 
(i.e. K — 0) and employed the Ito-discretization scheme with fixed time step dt, see Sec. MI B II In all simulation we 
have used an ensemble size of N — 1000 particles. A characteristic unit system was fixed by setting m — c = v = 1. 
Formally, this corresponds to considering re-scaled dimensionless quantities, such as p 1 = p 1 / (mc), x 1 = x t v/c, t = tv, 
V % = v i /c etc.. The simulation time-step was always chosen as dt = O.OOli^ -1 , and the Gaussian random variables 
y l {t) were generated with the pseudo-random number generator of Mathematica [44|. 



A. Distribution functions 



In the simulations we have numerically measured the stationary cumulative distribution function F of the absolute 
velocity values v = \JviV l in the laboratory frame So. Given, e.g., a spherically symmetric probability density 
4>{v) = (f>(v) with normalization 

1 = J d 3 v <f>(v) = 4tt^ dv v 2 j)(v), (50) 

the respective cumulative distribution function is defined by 

F(v) = 47T / du u 2 4>(u). (51) 
Jo 

In order to obtain F(v) from numerical simulations, one simply measures the relative fraction of particles with absolute 
velocities in the interval [0, v). Figure^depicts the numerically determined stationary distribution functions, taken at 
time t = lOOf -1 and also the corresponding analytical curves -Fi/sf/hk(^)- The latter were obtained by numerically 
integrating the formula (|51|l . using the three different stationary density functions </'i/sf/hk(' u ) = 0i/sf/hk( u ) found 
in Sec. HTH 

As one can see in diagram ^ (a), for low temperature values, corresponding to x 2> 1, the three stationary 
distribution functions approach each other, since they all converge to the nonrelativistic Maxwell distribution in 
the limit x — > oo. For high temperatures, corresponding to x 5~ 1> the stationary solutions exhibit significant 
quantitative differences, cf. Fig. ^(b). Since our simulations are based on an Ito-discretization scheme, the numerical 
data points agree best with the Ito solution (solid line). Similar to the (l+l)-dimensional case 0, the quality of the 
fit remains satisfactory over several orders of magnitudes of the parameter x- One may therefore conclude that the 
numerical simulations of Langevin equations provide a useful tool for studying relativistic Brownian motions in (1+3) 
dimensions. It should, however, be stressed again, that the appropriate choice of the discretization rule is particularly 
important with regard to potential applications to physical situations. 

Furthermore, Fig. (2J illustrates the influence of the number of spatial dimensions on the occurrence of relativistic 
effects. In the two diagrams we depicted the ratio Fhr(v) / Fm(v) for three different values of the characteristic pa- 
rameter Xi with .Fhk/mM denoting the cumulative velocity distribution function of the relativistic and nonrelativistic 
Maxwell distribution, respectively. Figure[5](a) corresponds to the case of (l+3)-dimensional free Brownian motions, 
while Fig. [3(b) refers to the (l+l)-dimensional case, discussed in paper I. The comparison of the two diagrams reveals 
that relativistic effects become significantly enhanced at lower temperatures (or larger values of Xi respectively), if 
the Brownian particle moves in (l+3)-dimensions. 
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0.2 0.4 0.6 0.8 1 0.2 0.4 0.6 0.8 1 



V [c] V [c] 

FIG. 1: These diagrams depict a comparison among the numerical Ito-prrescription and the three analytical results for the 
corresponding stationary cumulative distribution function F(v) in the laboratory frame So- (a) In the nonrelativistic limit 
X 3> 1 the stationary solutions of the three different FPE are nearly indistinguishable, (b) In the strong relativistic limit 
case x < 1j however, the stationary solutions exhibit deviations from each other. Because our simulations are based on an 
Ito-discretization scheme, the numerical data points do agree with the Ito solution (solid line) . 
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0.1 



0.01 



(a) (1+3)-dimensional case 

0.2 0.4 0.6 
v [c] 



0.8 



10 



0.1 



0.01 



X=10.0 
X= 1.0 
X= 0.1 



* * * X X * 
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(b) (1+1)-dimensional case 

0.2 0.4 0.6 
v [c] 



0.8 



FIG. 2: The diagrams show the ratio Fhk(v)/Fm(v) for three different values of the characteristic parameter \ ~ me? /{k^T), 
where Fm and Fhk denote the cumulative distribution functions of absolute velocity for the nonrelativistic and relativistic 
Maxwell distributions, respectively. Note that the deviation between the nonrelativistic and the corresponding relativistic 
distribution functions are more pronounced (at fixed value of x) f° r the (l+3)-dimensional situation, i.e. in Fig. (a). Therefore, 
compared with the (1+1) -dimensional case in Fig. (b), relativistic effects become detectable best at finite, lower temperatures, 
if the Brownian dynamics proceeds in (l+3)-dimensions. 



B. Mean square displacement 



Next, the spatial mean square displacement of the free relativistic Brownian motion is investigated. Since this 
quantity is easily accessible in experiments, it has played an important role in the verification of the nonrelativistic 
theory. 

As before, we consider an ensemble of N independent, free Brownian particles with coordinates x^(t) in S and 
initial conditions x^(0) = 0, W(j)(0) = for i = 1, 2, . . . , N. The position mean value is obtained as 

1 N 

= *(*)(*)' (52) 

and the related second moment is given by 

^(*) = ^ I>w(*)] 2 ■ ( 53 ) 
»=i 
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The empirical mean square displacement can then be defined as follows 

a 2 (t)=^(t)-[x(t)} 2 . (54) 
Important results of the nonrelativistic theory of the three-dimensional Brownian motion read 

lim x(t) -> 0, (55a) 

t — > + oo 

li m f_W _ 3-2D x , (55b) 

where the constant 

^=m: = _5_ (56) 

is the nonrelativistic coefficient of diffusion in coordinate space [not to be confused with noise strength D = kT / {mv)\. 

It is therefore interesting to consider the asymptotic behavior of the quantity a 2 (t)/t for relativistic Brownian 
motions. In Fig. [21 (a) one can see the corresponding numerical results for different values of x? evaluated on the 
basis of the Ito-scheme from Sec. IIII B II As evident from this diagram, for each value of Xi t ne quantity a 2 (t)/t 
converges to a constant value. This means that (with respect to the laboratory frame So) the asymptotic mean 
square displacement of the free relativistic Brownian motions is again normal, i.e. it increases linearly with t. For 
completeness, we mention that according to our simulations the asymptotic relation l|55a|) holds in the relativistic 
case, too. 




FIG. 3: (a) Mean square displacement, divided by elapsed time t, as numerically calculated for different x~ va ^ nes i n the 
laboratory frame So (rest frame of the heat bath). As evident from this diagram, for the relativistic Brownian motion the 
related asymptotic mean square displacement grows linearly with t. (b) The coordinate space diffusion constant -Dfo(x) was 
numerically determined at time t = 50is~ 1 . The solid line corresponds to the empirical fitting formula D x (x) = c 2 f _1 (x + 6) _1 , 
which reduces to the classical nonrelativistic result D x ~ c 2 /(ux) = kT/(mv) for \ ^ 6. 



In spite of these similarities between nonrelativistic and relativistic theory, an essential difference consists in the 
explicit temperature dependence of the limit value 6D X . As illustrated in Fig. 0(b), the numerical limit values 6Z?f , 
measured at time t = 50^ _1 , are reasonably well fitted by the formula 

D x = — r— ^— pr, (57) 

Kx + 6) 

which reduces to the nonrelativistic result (|56J) in the limit case ^ > 6 (low-temperature limit case). It remains as 
an open problem for the future to find an analytic expression for the relativistic diffusion constant D x . Note also 
that, with the position being the integral over the velocity degree of freedom, the different discretization rules do not 
impact the asymptotic (long time) result for the position diffusion coefficient. 



SUMMARY 



The challenge of this work has been to extend our previous work [23 on (l+l)-dimensional relativistic Brownian 
motions to the (l+3)-dimensional case. To this end, we have introduced in Sec. [H]a (l+3)-dimensional relativistic 
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generalization of the nonrelativistic Langevin equations (LE). Analogous to the nonrelativistic Ornstein-Uhlenbeck 
theory of Brownian motion [3I l24l I25L l45j . it is implicitly assumed that the heat bath, (which causes the stochastic 
motions of the particle) can be regarded as an isotropic, homogeneous fluid. Based on this assumption, the relativistic 
equations of motions are constructed such that they reduce to the well-known nonrelativistic LE in the limit case 
c — > 00. 

In our relativistic version of the LE, the viscous friction between Brownian particle and heat bath is modelled by a 
friction tensor v a p, exhibiting the same formal structure as the pressure tensor of a perfect fluid 38]. In particular, 
this means that the friction tensor is uniquely determined by the value of the (scalar) viscous friction coefficient u, 
as measured in the instantaneous rest frame of the particle. Similarly, the amplitude of the stochastic force is also 
governed by a single parameter D, specifying the Gaussian fluctuations of the heat bath, as seen in the instantaneous 
rest frame of the particle. 

In Sec. IIIBI we have rewritten the relativistic LE in laboratory coordinates, corresponding to a specific class 
of Lorentz frames, in which the heat bath is assumed to be at rest at all times. Further, it was shown that the 
relativistic equations can be recast such that they contain ordinary 'multiplicative' Gaussian white noise. Analogous to 
nonrelativistic stochastic processes with 'multiplicative' noise, this leads to an ambiguity regarding the interpretation 
of the stochastic differential equation' i.e., different discretization rules yield different Fokker-Planck equations (Sec. 
IIIIIl . Similar to the previous paper [23], we concentrated here on the three most popular discretization schemes, 
corresponding to the pre-point rule proposed by Ito PH|23| . the mid-point rule by Stratonovich and Fisk |23 . l29ll30ll3l| . 
and the post-point rule of Klimontovich and Hanggi (33, OS, yj, |35| . It was then shown in Sec. IIII Gl that only the 
latter prescription, i.e. the post-point discretization scheme, yields a Fokker-Planck equation, whose statio nary 
solution coincides with the relativistic Maxwell-Boltzmann distribution, as e.g. known from the work of Juttner |36| . 
Schay @ and de Groot et al. |22j . 

In Sec. II VI we presented several numerical results (based on the Ito scheme), including numerically obtained velocity 
distribution functions and furthermore, the mean square displacement of free Brownian particles. According to our 
findings, the relativistic mean square displacement grows linearly with the laboratory coordinate time; compared 
with nonrelativistic diffusions, however, the temperature dependence of the spatial diffusion constant becomes more 
intricate. 

Finally, we would like to emphasize that, so far, our approach of constructing a relativistic Brownian motion 
dynamics is merely based on the condition that - given the (a priori prescribed) stochastic force of the heat bath 
- the relativistic LE has to converge to the well-known nonrelativistic dynamical equation, if the Brownian particle 
moves sufficiently slow. In particular, there presently remains the open (and seemingly very difficult) problem of 
how to tackle in a relativistically consistent manner the dynamics of all those particles that constitute the heat bath 
(including their coupling to the relativistic Brownian particle). 
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APPENDIX A: STATIONARY SOLUTIONS OF RELATIVISTIC FOKKER-PLANCK EQUATIONS 



We seek the stationary solutions f(p) of the FPE l|29|) . which lead to vanishing currents 

Ji/sf/hk(p) = 0- 

We use the following Ansatz 

f(p) = C7~ Q exp(-x7), 

where C > is a normalization constant and 

, PiP X 

7= U+ —22 



(Al) 
(A2) 

(A3) 



The parameters a and x have to be determined from the condition (|A1|) . Therefore we have the following partial 
derivatives 



dl _ Pj 



dpi "fm 2 c 2 ' 
9f _ ~Pj 



dpi 'ymPc 2 \ 7 



+ X f 



and furthermore the divergence 



dpi K ' dpi 



5 lJ 



p % p? \ 1 



m 2 c 2 J 7 



3p l 



"fm 2 c 2 



(A4a) 
(A4b) 

(A4c) 



1. Ito current 



For the Ito current j I from 1|33[) . the condition IjAlfl yields 

= -ilip) 

- ¥f + Df^(A-if + D(A-if§l, 



P l f 



"fm 2 c 2 
3D 



m 2 c 2 J "f 2 m 2 c 2 V 7 



' - + X ) f 



"fm 2 c 2 m^c* \7 



D fa 



which is fulfilled for a — 3 and x — vm 2 c 2 /D. 



2. Stratonovich-Fisk current 

For the Stratonovich-Fisk current j SF from l|38|l . the condition (|A1(I becomes 

= -j§ F (p) 

where has been given in (|25[1 . A lengthy, though straightforward calculation shows that 



Inserting this into i|A6(l . one finds 



3p 4 



= vp l f + Df- ., ., 
Aymrcr 



D S ij 



p l p> 



-Pj ( a 



9 ? I ? 9 2. 



+ X f 



p l f 



3D 



D fa 



2~/m 2 c 2 m 2 c 2 \"f 
which is fulfilled for a — 3/2 and x — vm 2 c? /D. 



3. HK-current 

For the HK-current j K from l|41|l , the condition i|Al|l yields 

= -jkM 

= vp if + D{A- i y^ 



p l p 3 



vp*f + D[S«+ ., , , ., ., , 

m c J j^m^c^ \7 



-Pj 



- + x)f 



P l f 



which is fulfilled for a = and x = vmrcr jD 
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